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Abstract. We have compared the cosmic microwave background (CMB) temperature anisotropy maps made from one-year 
time ordered data (TOD) streams that simulated observations of the originally planned 100 GHz Planck Low Frequency 
Instrument (LFI). The maps were made with three different codes. Two of these, ROMA and MapCUMBA, were implemen- 
[ tations of maximum-likelihood (ML) map-making, whereas the third was an implementation of the destriping algorithm. The 

■~>^ ■ purpose of this paper is to compare these two methods, ML and destriping, in terms of the maps they produce and the angular 

' power spectrum estimates derived from these maps. The difference in the maps produced by the two ML codes was found to 

be negligible. As expected, ML was found to produce maps with lower residual noise than destriping. In addition to residual 
^ ■ noise, the maps also contain an error which is due to the effect of subpixel structure in the signal on the map-making method. 

' This error is larger for ML than for destriping. If this error is not corrected a bias will be introduced in the power spectrum 

, estimates. This study is related to Planck activities. 
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1. Introduction Frequency Instrument) channel of the Planck satellite. This 



simulated data contains contributions from CMB and fore- 
ground emissions, as well as from the instrumental noise. The 
noise was assumed Gaussian distributed and stationary over the 
entire mission, containing a 1// and a white noise component. 
We considered temperature anisotropics only; no polarisation. 



Map-making from an observed time ordered data (TOD) stream 
is an important step in the data processing pipeline of a cosmic 
microwave background (CMB) experiment. A number of map- 
making algorithms which, under the assumption of Gaussian 
distributed and stationary noise, aim at finding the optimal 

minimum variance map have been proposed (Wright 19961 Output maps were generated from the TODs using three 

Borrill .1999. Dore et al. .200 L Natoli et al. .2001) . The de- distinct map-making codes. The output maps and their angu- 

striping technique (Burigana et al. 119971 Delabrouille 119981 lar power spectra were critically compared. The map-making 

Maino et al.ll999 2002" Revenu et al.'2000al l2000bl Keihanen codes were ROMA (Roma Optimal Mapmaking Algorithm; 

et al. l2004> is a simpler map-making method. NatoU et al. 2001, de Gasperis et al 2005 ). MapCUMBA (orig- 

For this study we utihzed simulated one-year TOD streams, inally introduced by Dore et al. 1200 II current version based on 

that had been produced in the course of the work of the CTP the preconditioned conjugate gradient method) and destriping 

(Ce for Temperature and Polarisation) Working Group of the (Keihanen et al. ,2004j . All methods have been developed to 

Planck Consortia. These TODs represent the output of a sin- treat PLANCK-Uke data. 

gle detector of the originally planned 100 GHz LFI (Low a ■ ■ ■ • • i i l u j f 

f ^ ^ ^ A minimum variance map maximizes the likelihood func- 

Send offprint requests to: T. Poutanen, e-mail: tion involving the full noise covariance. The ROMA and 

torsti . poutanenOhelsinki . f i MapCUMBA algorithms aim at producing the minimum vari- 
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ance map. To accomplish this, the algorithms require knowl- 
edge on the characteristics of the instrument noise. Both iter- 
ative (Dore et al. 2001 1 and non-iterative (Natoli et al. 2002 1 
methods to estimate the noise properties directly from the data 
have been proposed. In this paper ROMA and MapCUMBA are 
referred to with a common name maximum likelihood (ML) 
map-making. 

Destriping does not employ the noise covariance matrix 
and does not aim at producing a minimum variance map in 
this sense. Thus it does not require prior knowledge on the 
characteristics of the instrument noise. This simplifies the al- 
gorithm considerably as compared to the ML map-making. 
In spite of this, destriping is able to provide an estimate of 
the low-frequency part of the instrument noise and to return 
a TOD where these noise components have been removed. 
Destriping can also be applied to estimate various systematic 
effects and drifts, remove them and return a cleaned TOD (see 
e.g. Mennella et al. l2002t . 

The aim of this paper is to compare these two methods, 
ML (ROMA and MapCUMBA) and destriping, in terms of the 
maps they produce and the angular power spectrum (Cf) esti- 
mates derived from these maps. 

This paper is organized as follows. In Sect. |2] we describe 
the map-making and power spectrum estimation methods ap- 
plied in this study. The simulated one-year TOD streams are 
introduced in Sect.|3l The output maps are compared in Sect.|3 
and the power spectrum estimates produced from the simulated 
TODs are examined in Sect. |5] The conclusions are given in 
Sect.|6l In AppendixlAlwe describe how the output map is split 
in the (wanted) binned noiseless map and in the (unwanted) re- 
construction error map. These quantities were considered in the 
comparison of the output maps in Sect. |3 In Appendix B we 
discuss some details about how the C( of the maps are related 
to the input C( used to generate the TODs. 

2. Methods 

Let us denote by a column vector y the samples of the observed 
TOD. The length of y is A^,, the number of samples in the to- 
tal mission. In the map-making problem we assume that the 
signal samples are scanned from a pixelized temperature map 
(m). (This assumption of course represents an approximation 
to reality, and thus contributes to error in the final maps.) The 
length of the column vector m is A^pix, the number of pixels in 
the map. 

The scanning is implemented by a pointing matrix P. The 
size of the pointing matrix is [P] = {Ni,Npix)- Each row con- 
tains zeros except for one element with value one, indicating 
the pixel at which the detector beam centre was pointing when 
the sample was measured. The map-making methods discussed 
here utilize pointing information only to the accuracy given by 
the pixel size. If the instrument beam response is spherically 
symmetric (and as we are neglecting polarisation), information 
about the rotation angle of the detector around the line of sight 
is not needed, and such a simple pointing matrix is sufficient. 
The pixel temperature of the output map will then represent a 
convolution of the sky and the beam. We have used in this study 
simulated TODs corresponding to both spherically symmetric 



and elliptic beams; but a study of the effects of beam shape is 
postponed to a future detailed study by the CTP group. 

For an asymmetric beam every pixel will be smoothed with 
a different response that is the mean of the beam orientations 
of the observations falling in that pixel. Recently a deconvo- 
lution map-making algorithm was introduced that can provide 
an output map where the smoothing of the instrument beam 
(symmetric or asymmetric) is deconvolved leading to a map 
that is an estimate of the true sky (Armitage & Wandelt [20041 . 
Deconvolution map-making is beyond the scope of this study. 

2.1. ROMA and MapCUMBA 

The output map (m) is solved by minimizing the log-Ukelihood 
formula (Natoli et sl. fIUUT\ 

=(y-Pm)^N-'(y-Pm). (1) 

Here N is the noise covariance matrix N = (nn^), where n is 
the instrument noise component of the TOD and (■) denotes 
expectation value. A set of linear equations is obtained for the 
output map 

P^N^'Pm = P^N^'y- (2) 

Eq. (|2j is a general result for the minimum variance map. 

Usually ML map-making assumes the noise to be stationary 
throughout the mission. It is further assumed that the elements 
of the covariance matrix (Nij) vanish when \i - j\ is larger than 
some A^,, and Njj <s: N,. This means that the correlation is signif- 
icant only across a number of samples that is a tiny fraction of 
the total length of the TOD. Thus the noise correlation matrix N 
can be approximated by a circulant matrix (Natoli et al. l2001> . 
Note that the matrix N"' is approximately circulant as well. 
The multiplication N^'y can be carried out more easily in the 
frequency domain where N"' is diagonal (Natoli et al. 12001 > . 
Both in the ROMA and in the MapCUMBA algorithms the out- 
put map is solved from Eq. (|2j with an iterative preconditioned 
conjugate gradient method. The iterations are repeated until the 
fractional difference has reached a low enough value (Natoli et 
al. 2001 1. This limit is typically on the order of 10"^. 

Due to the circulant matrix approximation each row of the 
matrix N"' contains the same element values with a different 
cyclic permutation. It is assumed that only the elements of a 
row with |/ - j\ < (Ng N,) have non-zero values. The rest 
of the elements are zero. The collection of the non-zero ele- 
ments (of a row) is called the noise filter. The lag of an element 
is the difference of its indices (i-j). The choice of the value A'^ 
is a significant decision for the quality of the output maps and 
for the computation time of the algorithm (Natoli et al. 200TJ. 

Natoli et al. ( 1200 U have shown how ML map-making can 
be applied when the instrument noise is only piece-wise sta- 
tionary. 

2.2. Destriping technique 

The destriping technique for map-making has been derived 
from the COBRAS/SAMBA Phase-A study (BersanelU et 
al. il996j . It has been implemented by several groups (Burigana 



Poutanen et al.: Comparison of map-making algorithms for CMB experiments 



3 



et al JTWn Delabrouille fT^ Maino et al JT^HM^ Revenu 
et al. I2000al l2000bl Keihanen et al. l2004l Efstathiouj^OOS I. 
The implementation studied in this paper (Keihanen et 
al. l2004j makes use of the fact that Planck is a spinning space- 
craft. Detector beams are drawing almost great circles on the 
sky. Each scan circle is observed several times before the spin 
axis is repointed. In order to reduce the level of instrumental 
noise, the signal can be averaged ("coadded") over these scan 
circles. In the following we call this averaged scan circle a ring. 
This coadding shortens the TOD by a factor given by the num- 
ber of circles between repointings. 

Janssen et al. ( 1996 1 has pointed out that the effect of the 
instrumental noise, in particular 1 // noise, on the ring can be 
approximated by a uniform offset or "baseline". The key prob- 
lem in destriping is to find the amplitudes of these baselines. 
The destriping technique uses the redundancy of the observing 
strategy by considering the intersections (crossing points) be- 
tween the scan circles to obtain these amplitudes. A crossing 
point is defined as two samples from different rings falling on 
the same pixel. According to this definition, two closely spaced 
rings may have a crossing point (or a sequence of them, when 
the rings run parallel) without actually crossing each other, if 
the distance between the rings is smaller than the pixel size. 

The correlated noise component of the TOD is modelled as 
Hcon = Fa (Keihanen et al. 2004). Here vector a contains the 
amplitudes of the baselines and matrix F unfolds them into a 
TOD. That is. Fa is a piecewise constant TOD, where the con- 
stants are given by the elements of a. Once the amplitudes have 
been solved. Fa is subtracted from the original TOD to pro- 
duce a cleaned TOD. Finally the output map is binned from the 
cleaned TOD. This procedure is formally described by min- 
imizing the following likelihood function with respect to the 
output map m and the amplitudes a (Keihanen et al. 120041 

X^^{y-Pm- Fa)^(y - Pm - Fa)/(rl (3) 

It is assumed here that the variance (0-2) of the non-correlated 
component of the noise is constant throughout the TOD. The 
amplitudes (a) can be solved from the equation 

F^ZFa = F^Zy (4) 

and the output map is given by 

m = (P^P)-'P^(y - Fa). (5) 

In Eq. (lUl Z s I, - P(P^P)"'P^, where I, is a unit matrix 
with dimension A^,. The matrix Z is determined by the cross- 
ing points of the scan circles. When Z is acting on the TOD it 
subtracts from each sample the average of the samples hitting 
the same pixel. 

It turns out that the matrix F^ZF is singular and additional 
conditions are required before the amplitudes can be solved 
from Eq. @. To accomplish this, we require that the sum of 
the baseline amplitudes is zero. 

The destriping method, as given by Eqs. @ and ^ requires 
no knowledge of the noise power spectrum. (If it is known that 
o"^ varies between different parts of the TOD, this information 
can of course be incorporated to improve the accuracy of the 
method.) 



The dimension of the matrix F^ZF equals the number of 
the fitted baselines. A typical number of baselines for a one 
year TOD (number of rings) is of the order of several thou- 
sands. This is a less complex problem than solving the map in 
the ML method (cf. Eq. Q), because the number of pixels in 
the output map is typically between 10^ . . . 10**. 

Generalized approaches to destriping method have been 
implemented (but not used in the present study) which are able 
to fit different sets of base functions (in addition to the uniform 
baseline) and may better remove the contributions of different 
systematic effects from the TODs (Delabrouille 19981 Maino 
et al. EUDl Keihanen et al. lMHI 120051 . 

2.3. Estimation of the angular power spectrum 

We studied the CMB angular power spectrum estimates ob- 
tained from the output maps. For the power spectrum estima- 
tion we used the MASTER approach (Monte carlo Apodised 
Spherical Transform EstimatoR) described by Hivon et al. 
(2002). MASTER has been adapted to operate both with 
ROMA (Balbi et al. I2002> and with destriping (Poutanen et 
al. 2004j . 

The angular power spectrum of the CMB sky Cf is derived 
from the spherical harmonic expansion coefficients (ae„) of the 
sky 

1 ^ 

= 27TT 2 (6) 

m=-e 

The simulated TODs are generated from such a spectrum, sup- 
posed to represent the "real sky". We call it the input spec- 
trum for our simulation + map-making problem. The ae„ coef- 
ficients are a realisation of the underlying "theoretical" angular 
power spectrum C*. The ensemble mean (of many realisations) 
equals the theoretical spectrum, {Cf) - Cf. The angular power 
spectrum (pseudo spectrum) of an output map is denoted Q. 
We further define C^, which is the angular power spectrum of 
the binned noiseless map. That map is obtained by binning the 
samples of the noiseless signal-only TOD to map pixels. 

The ensemble mean of the CMB pseudo spectrum depends 
on the ensemble mean of the spectrum of the sky. Hivon et 
al. l2002l express this relation as 

<Q) = 2 Ma,F(,BlDl{Cf) + (Ne). (7) 

c 

The matrix M(e> is the mode coupling matrix (kernel matrix) 
determined by the applied sky cut (Hivon et al. I2002> . A 
symmetric instrument beam is assumed and its smoothing is 
modelled by B^. Pixelisation introduces additional smoothing, 
which is represented by the pixel window factor D^. The fil- 
ter function F( represents a possible distorting effect of map- 
making and the noise bias (N() the remaining noise. 

There are a number of complications in the relation be- 
tween the spectrum of the output map C( and that of the sky 
Cf, not fully captured by Eq. 0. These are related to the ex- 
perimental setup, and, in the case of simulated data, to imper- 
fections in how the simulation models the experiment. Some of 
these are discussed in Appendix [5] The purpose of this study 
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is to compare map-making algorithms. We want to isolate the 
map-making errors from these other effects. Thus we will not 
try to estimate the angular power spectrum of the sky (C") but 
instead we compare the spectra of the output maps to the spec- 
trum of the binned noiseless map (C^). The map-making meth- 
ods compared here are derived from assumptions (see begin- 
ning of Sect. |2ji which correspond to the binned noiseless map 
being equal to the true sky (or the covered part of it), and thus 
it is effectively the object the map-making methods are trying 
to estimate. 

The spectra of the output map and the binned noiseless map 
can be related by 

<Q) = FtiC^) + m. (8) 

Here F[ accounts for the map-making errors only, and we shall 
use this definition for the filter function. It can be determined 
by e.g. signal-only MC simulations. We expect that its devi- 
ation from one should be small. We have dropped the mode 
coupling matrix because the sky coverage of the output map 
and the binned noiseless map are identical. Likewise, the beam 
(symmetric or not) and pixel window have the same effect on 
both maps, and therefore do not appear in Eq. (|8|l. The esti- 
mate of the spectrum of the binned noiseless map (Cf) can be 
obtained by inverting the equation 

Ct = F(Cf + (Ne) . (9) 

The influence of the instrument noise is modelled by the 
noise bias term {Ne}- For destriping an analytic method has 
been proposed (Efstathiou 2005 1 that can provide an estimate 
for the noise bias. However, in this study we used Monte Carlo 
(MC) simulations to obtain an estimate for it. A number of 
noise only TODs were generated from the power spectral den- 
sity (PSD) of the instrument noise. Maps were made from these 
TODs and a mean of their pseudo spectra was derived. This 
mean is an estimate of the noise bias. 

3. Time ordered data 

The "observed" TOD streams used in this study were gener- 
ated by computer simulations that used the Level S software 
(Reinecke et al. '2005). The correspondence between the sam- 
ple sequence of the TOD and locations on the sky is determined 
by the scanning strategy. The Planck satellite will be placed in 
an orbit around the 2""* Lagrangian point (L2) of the Earth-Sun 
system (Dupac & Tauber,20^j. A satellite placed around L2 
will stay near the ecUptic plane and will follow the Earth when 
it is orbiting the Sun. 

The Planck satellite rotates around its spin axis and the an- 
gle between the spin axis and the optical axis of the telescope 
is 85°. While the satellite is spinning (at nominal rate 1 rpm) 
the beam draws nearly great circles in the sky. The satellite spin 
axis is repointed at one-hour intervals. The diff'erent scanning 
strategies considered for Planck (Dupac & Tauber l2005> dif- 
fer in what path these repointings follow on the sky. Between 
the repointings the spin axis remains fixed. We applied a "cy- 
cloidal" scanning strategy where the spin axis followed a cir- 
cular path around the anti-solar direction. The angle between 




:^.5 ^^^^^K^ ^^^^H 4.0 



Fig. 1. Number of hits per pixel for the scanning strat- 
egy applied in this study. The map is in the ecUp- 
tic coordinate system. The scale is logi()(nhit), where 
nhit is the number of hits in a pixel. (A version of 
the paper with a better-quality figure can be found at 
http://www.physics. helsinki.fi/~tfo jcosm/tfo-planck.html ) 



the spin axis and the Sun-Earth axis was 10°, and the spin 
axis completed a full circle around the Sun-Earth axis every 
6 months (while the Sun-Earth axis itself of course made a full 
circle along the ecliptic in one year). 

The simulated TODs were generated for the originally 
planned 100 GHz LFI detector number 9. The length of the 
TODs was 12 months. They consisted of 525 960 scanning cir- 
cles with 6498 samples on each circle, corresponding to a sam- 
pling frequency of /s - 108.3 Hz. Since we assumed idealized 
satellite motion, where the 60 scan circles between repointings 
fell exactly on each other, sample by sample, these circles could 
easily be averaged into a single ring. This coadding of the TOD 
was performed before destriping was applied, but it was not 
done for the ML codes. 

In this study we utilized the HEALPix' pixelisation 
scheme. Its pixel dimension is set by the A^side resolution pa- 
rameter A map of the full sky contains 12N^.^^^ pixels. 

The number of hits per pixel (A^side = 512) of the applied 
scanning strategy is shown in Fig.[n At this resolution the sky 
coverage was 100 % (every pixel was hit). 

The instrument noise was a sum of white and l/f noise. 
The PSD of the noise was 

= + (/>/mi„). (10) 

Here /k is the knee frequency where the spectral powers of the 
l/f and white noise are equal. The nominal white noise stan- 
dard deviation (std) per integration time (f - l//s) is cr and 
/min is the minimum frequency below which the noise spec- 
trum becomes flat. The values of the noise parameters are given 
in Table [2 They represented a realistic expected noise perfor- 
mance of the instrument hardware. We used the stochastic dif- 
ferential equation (SDE) algorithm (from Level S) to gener- 
ate the TODs of the instrumental noise. Perfect knowledge of 

' http://healpix.jpl.nasa.gov 
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the noise parameter values was assumed both in the ML map- 
making and in the power spectrum estimation. 

Table 1. Simulation parameters used in the TOD generation. Two 
TODs (signal-only and signal+noise) were generated for all 4 sim- 
ulation cases leading to 8 TODs in total. The TODs were generated 
using Level S software (Reinecke et al.'2005'(. The TODs are sums of 
CMB (C), foreground (F), and/or noise (N) as indicated in the table. 
The CMB and foreground TODs were made using totalconvolver and 
interpolation algorithms (cases 1, 2 and 4) or they were scanned from 
a high resolution map (case 3). This is also indicated in the table. 

Parameters common to all cases 



Detector 

Number of detectors 
Mission time 
Scanning (a) 
Noise 



Parameters for Case 1 



<r(b) 

fm'm 



LFI 100 GHz 
1 

12 months 
Cycloidal 

3957.26 ^iK 
10-^ Hz 
108.3 Hz 



Case 2 Case 3 



Case 4 



TOD (S only) C+F (c) C+F (c) C (d) C (c) 

TOD(S-l-N) C+F+N C+F+N C-l-N C+N 

Beam (e,f) Symmetric Elliptic Symmetric Elliptic 

Noise - /k 0.03 Hz 0.03 Hz 0.1 Hz 0.1 Hz 

(a) 10° amplitude and 6 months period. 

(b) White noise std of the detector TOD (in antenna temp, units). 

(c) Signal made using totalconvolver and interpolation. 

(d) Signal scanned from a high resolution map. 

(e) Symmetric Gaussian beam: 

Full width half maximum (FWHM) = 10.6551 arcmin. 

(f) Elliptic Gaussian beam: 



FWHM,, 



11.8652 arcmin, FWHM,, 



: 9.5684 arcmin. 



We used 8 distinct TODs for this study. They were split in 
four cases (2 TODs in a case). The simulation parameters for 
the TODs are shown in Tabled The differences between the 
cases were the knee frequencies of the instrument noise, the 
telescope beams, whether the foreground^ was included, and 
the ways how the TODs were generated. The sky contained 
CMB and foreground emissions in cases 1 and 2 but only CMB 
in cases 3 and 4. Four TODs contained signal-fnoise and four 
TODs had signal only. The TODs were convolved with either a 
symmetric or an elliptic instrument beam (see Table^. 

The CMB signal was derived from a set of ac„, expansion 
coefficients that was a realisation from a theoretical CMB an- 
gular power spectrum C*. This was computed using the 
CMBFAST code^ (Seljak & Zaldarriaga 1996 1, and it cor- 
responds to a ACDM (cosmological constant + Cold Dark 
Matter) model. 

For the cases 1, 2 and 4, the expansion coefficients of the 
sky were convolved with the beam using the total convolu- 



^ Planck foreground template maps were applied in the TOD gen- 
eration. The templates are available for the Planck collaboration at 
j http://planck.mpa-garching.mpg.de 

^ http://www.cmbfast.org 



tion technique (Wandelt & Gorski EOOlt . The totalconvolver 
algorithm (part of Level S) outputs a discrete temperature field 
which is tabulated in an equally spaced three dimensional grid 
(two dimensions for the pointing of the beam centre and one 
dimension for the beam orientation). The TOD samples were 
interpolated from the tabulated temperature grid. For the case 

3, the signal TOD was scanned from a high resolution map gen- 
erated from the a(m that had been convolved with the beam. The 
map had /Vside - 1024 and it was generated with the SYNFAST 
code of the HEALPix package. 

4. Output maps 

We had 8 TODs available for map-making comprising 4 cases 
(a signal-i-noise TOD and a signal-only TOD for each case). 
Output maps were made from these TODs using three map- 
making codes, leading to 24 output maps in total. 

The 4 cases differed in a number of ways: in whether fore- 
ground was included, in noise, in beam shape, and in how 
TODs were generated. These differences represent both real ef- 
fects and imperfections (e.g., interpolation errors) in TOD gen- 
eration. They result in variations of map-making performance 
from case to case. However, these variations are not the ob- 
ject of this paper (The TODs differed in too many ways for a 
study of the effect of these differences to be conducted from just 
4 cases). Instead, the object is to compare map-making meth- 
ods. We have included all 4 cases available to us, mainly to see 
whether the results of comparison between methods stay con- 
sistent from case to case. They do. Some conclusions regarding 
the effect of foregrounds on the different map-making methods 
can also be drawn. 

All output maps were made with pixel resolution /Vside - 
512. The noise in the signal-t-noise TODs resembled the noise 
of a single LFI 100 GHz detector. Thus the maps are much 
noisier than they would be if they were made from a full set of 
24 detectors, by a factor of about V24. 

Before the ROMA and MapCUMBA output maps could be 
made the noise filters had to be produced (see Sect. 12. 1> . They 
were determined from the analytical model of the noise PSD 
(see Eq. ilQ\ ) and its known parameter values. The noise fil- 
ters were symmetric and had - 65537 elements at non- 
negative lags (lag > 0). The applied noise filters were iden- 
tical in ROMA and in MapCUMBA. The conjugate gradient 
iterations were continued until the fractional difference had de- 
creased to < 10"^. 

For destriping the amount of TOD was reduced by averag- 
ing 60 scan circles between the repointings. The baseline am- 
plitudes were solved exactly (no iterations) from Eq. (0} using 
as an additional condition that their sum is zero. 

The output maps of ML map-making and destriping can 
be divided into a binned noiseless map and a reconstruction 
error map (Tegmark . 1997a) . The binned noiseless map is the 
signal-only TOD binned to map pixels. The reconstruction er- 
ror map (e) is the (unwanted) deviation of the output map from 
the binned noiseless map. A goal of map-making is to minimize 
this error. The reconstruction error map can be further split into 
the signal (e^) and noise (s„) components: e - e,, + e„. These 
are discussed in detail in AppendixlAl As shown there, the sig- 
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Fig. 2. The ROMA output map for case 1 . The units are antenna 
Kelvins at 100 GHz. The maps for MapCUMBA and destrip- 
ing look similar This map is an output from one detector For 
this plot the map resolution was degraded to A^side = 256. The 
monopole was removed from the map. 





Fig. 3. The ROMA reconstruction error map for case 1. The 
units are antenna Kelvins. The maps for MapCUMBA and de- 
striping look similar. The map is an output from one detector. 
For this plot the map resolution was degraded to A^side - 256. 
(A version of the paper with a better-quality figure can be found 
at http://www.physics.helsinki.fi/~tfoxosm/tfo_plan ck.htmli ) 



nal component arises frompixelisation noise (Dore et al.'^OOTi. 
The noise component (e„) is the output map from the noise part 
of the TOD. 

Note that the map-making methods are linear (up to numer- 
ical accuracy effects), so that these components can be obtained 
separately and studied independently, by applying the codes to 
the noise and signal parts of the simulated TODs (see e.g. de 



Gasperis et al. l2005t . 

Our maps express the temperature fluctuations in the an- 
tenna temperature. The ratio of the thermodynamic tempera- 
ture fluctuation to the antenna temperature fluctuation is (e ' - 
l)^/x^e'", where x = hv/kTo, h is the Planck constant, v is the 
frequency, k is the Boltzmann constant and To - 2.725 K is 
the CMB temperature. For this study the ratio is 1.287 (v = 
100 GHz). 




Fig. 4. a) Difference map between the output maps of ROMA 
and destriping. b) Difference map between the output maps 
of ROMA and MapCUMBA. The pixel minimum, maximum 
and std values (from the A^side =512 maps) for map a) are (- 
109.7, 61.1, 10.0914) yuK and for map b) (-27.4, 38.7, 0.1531) 
yuK. The units are antenna Kelvins. All maps are for case 
1. The corresponding difference map between MapCUMBA 
and destriping looks similar to map a). For the plots the 
map resolution was degraded to N^i^^ — 256. (A version 
of the paper with a better-quality figure can be found at 
|http://www.physics.helsirLki.fi/~tfojcosm/tfo-planck.html| ) 



The output maps and the reconstruction error maps from 
all three codes look similar. Figs. 12 and |3 show them for case 
1. To see the difference, one needs to calculate the difference 
map between the output maps of the different codes (Fig. |4}. 
The prominent feature in the ROMA - MapCUMBA difference 
map (see Fig.|^) falls on top of the last repointing period of the 
scan. In MapCUMBA the TOD is extended by copying sam- 
ples from its beginning to the end and vice versa. In ROMA the 
TOD is extended by padding zeros to the end. The TOD ex- 
tension is required to make its length suitable for convolution 
with the noise filter (number of TOD samples to be an inte- 
ger multiple of 2(A'^ - 1)), but the added samples are not pro- 
jected on the output map. The feature appearing in the ROMA 
- MapCUMBA difference map reflects the different treatments 
of the end of the TOD. 

Reconstruction error maps (s) were made by subtracting 
the binned noiseless map from the signal-nnoise output maps. 
Minimum, maximum and std values of the pixel temperatures 



Poutanen et al.: Comparison of map-making algorithms for CMB experiments 



7 



Table 2. The std, minimum and maximum values (all in /jK at 
100 GHz antenna scale) of the pixel temperatures of the reconstruc- 
tion error maps (e = e,, + e„). The numbers were calculated from 
A'side = 512 maps. They can be compared to the std of a white noise 
map, calculated as ^[{TJrCja, where ( 1 jn) is the mean (taken over the 
hit pixels) of the inverse of the number of hits in a pixel. The val- 
ues are 137.219 ^JK for cases 1 and 2 and 137.400 /jK for cases 3 
and 4. The values are different between the cases due to small differ- 
ences in the cycloidal scannings. This white noise std represents the 
level below which one cannot get. The difference between ROMA 
and MapCUMBA is negligible. We give the numbers in the table with 
many digits just to show at what level this difference is. The numbers 
in this table and the next one should be compared in the horizontal 
direction, to see the difference between the codes. The differences in 
the vertical direction reflect the eff'ect of several differences in how the 
TODs were generated, but these effects are not the object of this paper. 



std for e 








(min, max) 


ROMA 


MapCUMBA 


Destriping 


Case 1 


138.2497 


138.2496 


138.455 




(-825.3, 876.0) 


(-825.3, 876.0) 


(-822.5, 891.1) 


Case 2 


138.2475 


138.2474 


138.454 




(-825.0, 876.7) 


(-825.0, 876.7) 


(-822.5, 891.2) 


Case 3 


138.9109 


138.9110 


139.397 




(-838.4, 938.2) 


(-838.1,938.1) 


(-851.2, 952.3) 


Case 4 


138.9114 


138.9114 


139.398 




(-838.1,937.9) 


(-838.1,937.9) 


(-851.2, 952.2) 



Table 3. Same as TableElbut now the values are for the signal com- 
ponent of the reconstruction error map (Bp). The units are antenna /jK 
(at 100 GHz). 



std for Bp 








(min, max) 


ROMA 


MapCUMBA 


Destriping 


Case 1 


0.8794 


0.8796 


0.281 




(-100.7,53.8) 


(-100.7, 53.8) 


(-3.0, 2.1) 


Case 2 


0.6140 


0.6142 


0.252 




(-62.2, 37.8) 


(-62.2, 37.8) 


(-2.3, 1.8) 


Case 3 


0.3130 


0.3130 


0.120 




(-2.0, 1.9) 


(-2.0, 1.9) 


(-0.6, 0.7) 


Case 4 


0.4064 


0.4064 


0.169 




(-2.8,3.1) 


(-2.8,3.1) 


(-0.9, 0.8) 



of the reconstruction error maps are given in Table |2l In terms 
of the map variance ML map-making is slightly better than de- 
striping. Table|2lalso shows that the ML std is higher than the 
std of the white noise indicating that excess noise remains in the 
map. The cases 3 and 4 have higher map std's than the cases 
1 and 2. This is mainly caused by the higher knee frequency 
of the instrument noise in the cases 3 and 4. The map vari- 
ances of ROMA and MapCUMBA are practically equal show- 
ing that their performances are similar. The results shown in 
Table|2]represent the performance of a single LFI detector For 
frequency maps made from the observations of multiple LFI 
detectors the noise levels would be correspondingly lower. For 
comparison it can be noted that the noise level of a one-year 



W-band (94 GHz) frequency map (A^^ide =512) of the WMAP^ 
experiment (from a total of eight W-band detectors making up 
4 differencing assemblies) is slightly higher (~142 yuK) than 
the noise levels of Table |2l which represents just a single LFI 
detector. 

Pixel statistics for the signal components of the reconstruc- 
tion error maps (e^,) are given in Tabled These maps were pro- 
duced by subtracting the binned noiseless map from the noise- 
less signal-only output maps. In ML map-making the TOD was 
convolved with the noise filter As an example, the signal com- 
ponents of the reconstruction error maps for case 1 are shown 
in Fig. |5] The maps contained CMB and foreground signals. 
The error magnitudes are larger for the ML map-making than 
for the destriping. Additionally, the largest errors in the ML 
map tend to locate where the foreground signal (actually, its 
gradient) is strongest. For destriping this error appears as (er- 
roneous) baselines corresponding to offsetting an entire ring. 
Therefore a similar correlation between the errors and the fore- 
ground signal is not visible in the destriped maps, although we 
expect that the largest eiTors still originate from where the fore- 
ground signal is strong. 

The methods that were used to generate the signal parts 
of the TODs were not uniform (total convolution vs. scanning 
of a high resolution map, see Sect. O. That complicates the 
comparison of the results of Table |3l between the beams. This 
comparison was not attempted in this study. However, the re- 
sults can well be compared (at a given beam) between different 
map-making algorithms which was the main purpose of this 
study. 

Table |3] shows that the signal-only ML maps deviate more 
from the binned noiseless map than what happens in destriping. 
Foreground increases this error and large errors may occur in 
some pixels of the ML maps. The std of Ep is small compared 
to the std of the overall reconstruction error map (see Table |2j 
which is an indication that the total error is dominated by the 
instrument noise. 

The source of Bp is the pixelisation noise (as shown in 
Appendix IaI. In ML map-making the pixelisation noise spec- 
trum up to the knee frequency of the instrument noise con- 
tributes to Ep. In destriping, however, only a lower frequency 
part contributes leading to a smaller Ep. This is addressed in 
more detail in Appendix|X] The galactic foreground signal has 
a stronger spatial variation than CMB. This results in higher 
pixelisation noise, which explains the higher Ep magnitudes in 
the cases 1 and 2 than in the cases 3 and 4 (see Table|3}- 

We calculated the angular power spectra of the recon- 
struction eiTor maps (e) (Fig. |6j. The spectra of ROMA and 
MapCUMBA were very similar and could not be distinguished 
in the plot. The ratio of the power spectra between destriping 
and ROMA is shown in Fig.0 It seems that in destriping e has 
higher power in most multipoles. 

To examine this further, 100 MC noise-only TODs were 
produced from the known PSD of the instrument noise (see Eq. 
MQV ). The MC noise TODs were generated by the SDK method 
using a different seed value for every realisation. To have a rea- 
sonable calculation time for this MC study, multiple 35 hour 

http://map.gsfc.nasa.gov 
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-1.0 ^^^^B. 1.0 ^/.K 



Fig. 5. The signal components of the reconstruction error maps 
for ROMA (a) and for destriping (b). These are difference 
maps between the signal-only output maps and the binned 
noiseless maps. The output maps contained CMB and fore- 
ground. The units are antenna Kelvins. All maps are for 
case 1. The corresponding map for MapCUMBA looks sim- 
ilar to map a). The map resolution of the plots is A^side = 
512. Note the different colour scales of the maps. (A ver- 
sion of the paper with a better-quality figure can be found at 
^http://www.physics. helsinki.fi/~tfoxosm/tfo_planck.html ) 

chunks of noise were generated simultaneously in parallel pro- 
cessing and the chunks were glued one after another at the end. 
This leads to an MC noise TOD with no correlation between 
the chunks. Because the correlation of the noise in the observed 
TODs is weak at > 35 hour lags^, the zero correlation beyond 
35 hours in the MC noise is expected to cause an insignificant 
error in the noise bias estimates. 

Output maps for ROMA and destriping were made from 
the MC noise TODs and their angular spectra were derived. 
The mean spectra are shown in Fig.|8] It shows that the mean 
angular power is higher in destriping at all scales. The lower 
power at some t (in Fig. Q seems to be just due to random 
variation. 

The map-making codes were run on an IBM SP RS/6000 
computer with a cluster of Power3 processors running at a 
clock speed of 375 MHz. ROMA and MapCUMBA codes were 

^ R(t)/R(0) k 4 X 10"\ where R(t) is the autocorrelation of the 
instrument noise (an inverse Fourier transform of the noise PSD Eq. 
<10t ) and T = 35 hours. 



run parallel in multiple processors (number of processors was 
typically between 192 and 256) and it took ~ 10 min to produce 
an ML output map. In destriping an output map was produced 
in ~7 min in a single processor job. In MC studies this time can 
be reduced to ~4 min by inverting the matrix (see Eq. @) once 
and using the inverse in the subsequent runs. 

Note that a part of this large difference in computation cost 
is due to destriping being applied to a coadded TOD, which was 
a factor 60 shorter than a full TOD. This coadding was possible 
without eiTor, because we assumed idealized pointing. In real- 
ity, the pointings of the different circles of the same ring do not 
fall exactly on top of each other. This means that either some 
additional error is introduced by the coadding, or that destrip- 
ing has to work with the actual pointings of the full TOD. The 
latter option increases the computational cost of destriping, but 
it will still be significantly less than for ML map-making. 




Fig. 6. The angular power spectra of the reconstruction error 
maps (e) for ROMA and destriping. The angular power spectra 
for ROMA and MapCUMBA were very similar and would be 
on top of each other in this plot. The curves are for case 4. 



5. Power spectrum estimates 

To see how the differences between destriping and ML map- 
making are reflected in the angular power spectrum estimates, 
we derived C[ estimates from the output maps of case 3, where 
the TOD contained CMB and noise. Since the output maps 
from ROMA and MapCUMBA were practically identical, we 
let ROMA represent ML map-making in this section. The CMB 
was convolved with the symmetric beam (case 3, see Table 
As defined in Sect. 12.31 our power spectrum estimates are 
estimates of the angular power spectrum of the binned noise- 
less map. The estimates were compared to the actual spectrum 
(C^) of the binned noiseless map. That map was binned from 
the noiseless TOD containing CMB only. The angular spectrum 
is shown in Fig. |9] which also shows the input spectrum 
C". The difference between and O" is due to a number of 
effects, which are discussed in Appendix IbI but they are not 
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Fig. 7 

curve. 



The ratio of the angular power spectra of Fig.|6l a) Full 
b) Zoom to low {. 



ROMA 
DESTRIPING 



10 




10 



10 



10 



Fig. 8. Mean angular power spectra of output maps from 100 
MC realisations of noise TODs. The realisations were derived 
from the known PSD of the instrument noise (see Eq. J10» . 
Because the MC TODs were noise-only, this plot is for cases 3 
and 4. 



relevant for our comparison between map-making methods, as 
they can only produce estimates of . 

The relation between the pseudo spectrum Ct (angular 
power spectrum obtained from the output map) and the power 
spectrum estimate (C^) was given in Eq. The estimate is 
obtained by inverting the equation. The estimate of the noise 
bias {Nc) was obtained from the MC simulations (see Fig.|8j- 
The value of F[ was initially set to one. The obtained power 
spectrum estimates, the spectrum of the binned noiseless map 
and the MC noise bias (for destriping) are shown in Fig. 1101 

For the quality of an angular power spectrum estimate its 
bias and covariance matrix are important figures of merit. For 
the covariance matrix we restricted our study to its diagonal 
elements, which represent the error bars of the estimate. 



200 400 600 800 1000 1200 1400 

/ 

Fig. 9. The angular power spectrum (C^) of the binned noise- 
less map (black curve). The map (A^side = 512) was binned 
from the simulated noiseless TOD containing CMB that was 
smoothed with the symmetric beam (case 3, see Table 0. The 
spectrum shown in the plot is deconvolved with the same 
symmetric beam response. For comparison, the input power 
spectrum CJ," of the CMB sky is shown as well (grey curve). 
The difference between and C™ is (in this case) mainly due 
to pixel window smoothing and is discussed in Appendix 151 

5.1. Bias 



We defined the estimation error ACc as the difference between 
the power spectrum estimate and the spectrum of the binned 
noiseless map: AC( - C^-Cf. The estimation error was binned 
by averaging Af multipoles to a bin 



= — y 



+ l)ACe 



2jt 



(11) 



We evaluated the binned errors ACt for ROMA and destriping 
from the spectra shown in Fig. [TO] To facilitate the comparison 
of the binned estimation errors, we normalized them by divid- 
ing them with an analytic approximation of their std, which was 
obtained from 



-y 



'+1) 



2n 



1/2 



where 



A/(2^+l)/,ky 



(12) 



(13) 



is the approximation for the std of the unbinned error AC( 
(Ef stathiou l2005> . In this formula Cf is a given signal (no cos- 
mic variance). Sky coverage fraction is /sky = 1 in our case. For 
all normalizations we used the {N(}mc obtained for destriping 
in place of {N(). 

The normalized errors are shown in Fig. ^2 Their std (from 
one noise realisation to another) should be ~1. If an angular 
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Fig. 10. a) Angular power spectrum estimate Cf for destriping 
(grey curve). The ROMA estimate would fall nearly on top of 
the destriping estimate and would not distinguish in this plot. 
The estimates were derived from a single sky realisation (case 
3). The output maps covered the full sky and contained CMB 
and instrument noise. Since just a single detector was consid- 
ered, the noise becomes dominant already at around £ ^ 350. 
Cf is the angular spectrum of the binned noiseless map and 
{Nt)MC is the MC noise bias (for destriping). Note that they 
have been deconvolved with a symmetric beam that was iden- 
tical to the instrument beam. The filter function F( was set to 
one for both estimates, b) Same as a) but the spectra have been 
£ binned to Ae - 25. The difference between the ROMA and 
destriping estimates is visible at ^ > 1000. 



power spectrum estimate has a non-zero bias the mean of the 
fluctuations of the normalized error will have a positive or neg- 
ative trend. In the case of zero bias the mean will be close to 
zero. We can note that the normalized errors are different for 
ROMA and destriping, the largest differences being at ^ < 800. 
The differences are, however, smaller than the std of the errors. 
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Fig. 11. The normalized power spectrum estimation errors for 
ROMA and destriping. The normalized error was obtained by 
dividing the binned estimation error (Eq. il l\ . Ad - 25) with 
the analytic approximation of its std (Eq. il2\ ). The std of 
destriping was used in all normalizations. The filter function 
had value 1.0 for both curves. The bin-to-bin fluctuations are 
mainly caused by the instrument noise. 
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Fig. 12. Filter functions F[ for ROMA and destriping. They 
were obtained from Eq. (|9} (with {N(} = 0) by dividing the 
pseudo spectra of the output maps of the noiseless CMB-only 
TOD with the spectrum of the binned noiseless map. The fil- 
ter function for MapCUMBA was nearly identical to the filter 
function of ROMA and would not distinguish from the ROMA 
filter function in this plot. 

We could expect that some bias could be introduced to our 
power spectrum estimates because we neglected (by setting Ff 
= 1 in Eq. (|9}) the error that the map-making causes to the 
CMB signal. This is a reflection of the signal component Sj, 
of the reconstruction error that we found in the map domain 
(see Sect.l^. To assess the level of the bias we estimated Fg for 
ROMA and for destriping. We carried out no MC simulations to 
estimate them, but we determined them from Eq. (|9j using the 
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Fig. 13. Same as Fig.[^but now the power spectrum estimates 
have been corrected with the estimates of the filter function 
(from Fig. ll2t . 



pseudo spectra C( from the output maps of the noiseless (CMB- 
only) TOD {{Ne} = 0). Because the values of Fc are based 
on one CMB realisation only, these results should be taken as 
indicative. (In fact, we now fully correct for the effect of Bp, 
since the filter function is derived from the same realisation to 
which it is applied. In reality, of course, the signal-only TOD 
will not be available, and the filter function should be evaluated 
as an expectation value. It will then remove only the bias due 
to Bp.) 

The obtained F( are shown in Fig. [21 For destriping there 
is essentially no filter function {F[ ~ 1, Poutanen et al. 2004 1. 
For ROMA the deviation from 1 is larger, showing its largest 
values at £ = 800 . . . 1000. If not corrected, the map-making 
errors cause a bias in the ML spectrum estimates whose max- 
imum value in this case would be ~0.6% of the magnitude of 
the CMB spectrum. 

We corrected our angular power spectrum estimates with 
the obtained Fg and reproduced the normalized estimation er- 
rors. The result is shown in Fig.^] When comparing to Fig. II II 
the improved match between the normalized errors of ROMA 
and destriping can be noted (especially at ^' < 600). The re- 
maining differences are mainly due to the differences in the 
noise of the output maps of these two algorithms. 

5.2. Error bars 

An error bar is defined here as the square root of the diag- 
onal element of the covariance matrix (ACfACf) (+lcr er- 
ror bar). The error bars can be derived either analytically 
(Tegmark '1997b Efstathiou EUO?! or by MC simulations 
(Hivon et al. 2002 Poutanen et al. 120048 . In this study we 
did not do signal-nnoise MC simulations to determine the er- 
ror bars, but used instead an approximation to compare ROMA 
and destriping. We used cr( from Eq. (I13> with a modification 



that takes into account the different filter functions for different 
map-making algorithms {F( assumed value 1.0 in Eq. il3i ) 



4 



<(AQ)2) 



(2^+1) A 



F^ 



(14) 



The spectrum represents here a given signal (no cosmic 
variance). 

Applying the spectrum from Fig. |5] the noise biases 
{{N[}mc) from Fig.|S]and the filter functions from Fig. 1121 the 
ratio of the error bars between destriping and ROMA was eval- 
uated. It is shown in Fig. (black curve). The error bars at 
( < 1000 are larger for destriping than for ROMA. The largest 
relative differences are ~5%. The main cause of the larger error 
bars is the higher level of noise in the output maps of destriping 
(see Table|2}- At high-f the larger map noise of the destriping 
is partly compensated by its larger filter function (see Fig.ll2> 
leading to error bars that have nearly the same magnitude as the 
ROMA error bars. 

As a second case we assumed that we want to estimate 
the angular spectrum of the underlying theoretical CMB 
(instead of its particular realisation as above). Cosmic vari- 
ance then increases the error bars (Scott et al. .1994> Hobson 
& Magueiio [T996t : 



<(AQ)2) 



2 FtCf + {Ne}MC 



V(2^+l)" 



(15) 



The ratio of these error bars is shown in Fig. as well (grey 
curve). Because noise dominates the error bars at high multi- 
poles, these are similar to the error bars without cosmic vari- 
ance. Due to the dominance of cosmic variance at low mul- 
tipoles, the magnitudes of the error bars of the two methods 
approach each other at low {. 

6. Conclusions 

We have presented a comparison of the maps produced by three 
different map-making codes and two map-making methods, de- 
striping and ML map-making. We also compared the angular 
power spectrum estimates obtained from destriping and ML 
maps. The maps and power spectra were derived from a set 
of one-year TOD streams that resembled the observations ex- 
pected from a single 100 GHz Planck LFI detector. 

In terms of the map variance the two ML codes, ROMA and 
MapCUMBA, produce nearly identical maps, with lower noise 
than destriping. This lower noise is an advantage for them and it 
facilitates smaller error bars for the power spectrum estimates. 
The difference is, however, rather small. 

ROMA and MapCUMBA require knowledge of the power 
spectrum of the instrument noise, whereas destriping does not. 
In a real experiment the noise spectrum (if required) needs to 
be estimated from the observed data. Some estimation error 
can be expected which may increase the noise in the ROMA 
and MapCUMBA maps. Thus differences in the noise perfor- 
mance between ROMA/MapCUMBA and destriping may be- 
come smaller in a real experiment. A perfectly known instru- 
ment noise spectrum was assumed in this study. 
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Fig. 14. An estimate for the ratio of the error bars between de- 
striping and ROMA. Black curve is for the estimation of a par- 
ticular CMB realisation of the sky (no cosmic variance, error 
bars from Eq. il4i ) and the grey curve is for the estimation 
of the underlying theoretical CMB spectrum (cosmic variance 
included, error bars from Eq. il5i ). 

The map-making methods caused errors (exhibited in the 
signal component of the reconstruction error map) in the signal 
part of the output maps. The origin of these errors is the sub- 
pixel structure of the signal (pixelisation noise). It was shown 
that these errors were smaller in destriping than in ROMA and 
MapCUMBA. It was further shown that, if a proper correction 
is not applied, these errors may show up as an extra bias in the 
power spectrum estimates. The extra bias would be larger for 
ROMA and MapCUMBA than for destriping. 

In terms of CPU resources destriping is less demanding. 
This is an advantage in e.g. MC simulations. 
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Appendix A: Reconstruction error map 

In this Appendix we discuss in more detail the reconstruction 
error map and how its signal and noise components arise. 



The output map of the ML map-making can be solved from 
Eq. (|2}. That equation is reproduced here 

P^N^'Pm = P^N^'y- (A.l) 

The noise covariance matrix N can be freely normalized 
by a constant without affecting the output map. Let us assume 
that each element of N has been divided by the variance (cr^) 
of the non-correlated noise component of the observed TOD 
(vector y). By replacing N ' with an identity I, - (I, -N~') and 
reaiTanging some of the terms one obtains for the output map 

m = (P^P)"'B"'P^N"'y, (A.2) 
where 

B = I,„ - P^(I, - N-')P(P^P)-'. (A.3) 

Here I, and I,„ are unit matrices with dimensions equal to the 
number of samples A^, of the TOD and the number of pixels 
A^pix in the map, respectively. We can apply a geometric series 
trick (Tegmark ..l997aj to prove the following identity 

- P^(I, - N-')P(P^P)-']-'P^ = 
= P^ [I, - (I, - ' )P(P^P)" ' P^] " ' . ( A.4) 

Applying this in Eqs. iA.2\ and (IA.3> and noting the definition 
of the matrix Z (see Sect. 12. 2> one obtains for the output map 

m = (P^P)-'P^[I, + (N- I,)Z]-'y- (A.5) 

By adding and subtracting y on the right hand side and carry- 
ing out some arithmetic manipulations the output map can be 
expressed in the following form 

m = (P^P)"'P^[y- A]. (A.6) 

Vector A is solved from the linear equation 

[Z + (N-I,)-']A = Zy. (A.7) 

Let the complex valued matrix H {[H] = (A^,, A^,)) be an in- 
verse DFT (Discrete Fourier Transform) operator that converts 
frequency domain vectors to time domain (to TOD domain): 
y = Hy, where y is the frequency domain counterpart of y. The 
inverse operator to H is its Hermitian conjugate H' . We can 
assume that the matrix H is normalized to HH^ - H'H - I,. 
A square matrix A ([A] - (N,,N,)) in the TOD domain can be 
converted to a matrix A in the frequency domain: A - H ' AH. 

After converting both sides of Eq. (j^T) into frequency do- 
main, the vector A (frequency domain counterpart of A) can be 
solved from the equation 

[H'ZH + (N- I,)"']A = H^Zy. (A.8) 

The solution for the output map becomes then 

m = (P^P)"'P^[y- HA]. (A.9) 

Comparing Eqs. (IA.8> and iA.9i to the corresponding Eqs. (0} 
and (|5} of destriping the similarity of the output map solutions 
between ML and destriping can be clearly seen. 
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The observed TOD (y) contains two components: signal 
s and instrument noise n (see Sect. O. The term Zy (see Eq. 
JA.8» can now be split in two components 

Zy = Zs + Zn. (A. 10) 

Writing out the first term on the right hand side we obtain 

Zs = s - P(P^Pr'P^s. (A. 11) 

Apart from the sign, Zs is the pixelisation noise introduced by 
Dore et al. (2001 1. Pixelisation noise represents error that is 
caused by the discretization of the sky into pixels. Following 
Dore et al. (2001 1 we define pixelisation noise p = -Zs. 

The split of Zy in two components means that A is spUt in 
two components as well: A = + A„ . They can be solved from 



[H'ZH + (N - I,)"']A,, = -H'p 
and 

[H^ZH + (N - I,)"']A„ = H^Zn. 



(A. 12) 



(A.13) 



In destriping the amplitudes of the base functions are split anal- 
ogously: a = Sp + a„. The first component is determined by the 
pixelisation noise and the second one by the instrument noise. 

Next we insert y = s -)- n into Eqs. (IA.6> and (IA.7> . We 
obtain for the output map 



m = (P^P)"'P^s - (P^P)"'P^HA,, + 
H-(P^P)"'P^[n-HA„]. 



(A. 14) 



Ideally, we would like the output map of the map-making algo- 
rithm to be equal to the first term on the right hand side. It is 
called binned noiseless map in this study. The rest of the terms 
bring error They are represented by the reconstruction error 
map (Tegmark ll997a> 



e = m - (P^P) 'P^s 



(A.15) 



The reconstruction error map is comprised of a signal compo- 
nent {Bp) and a noise component (e„). 



S — Sp + . 

Ep - -(P^P) ip^HA^ 

Sn = (P^P)-'P^[n - HA„] 

For destriping we can write analogously 

Bp = -(P^P)-ip^Fap 



e„ = (P^P)-ip^[n-Fa„]. 



(A. 16) 
(A.17) 
(A.18) 

(A. 19) 
(A.20) 



The signal and noise components of the reconstruction 
error map were used extensively when the map-making al- 
gorithms were compared in Sect. |4] We studied the mini- 
mum, maximum and std values of their pixel temperatures. 
Additionally, we produced their angular power spectra and 
compared them as well. 



A.1. Map-making errors and pixelisation noise 

The purpose of this section is to give a qualitative explanation 
to the fact that the magnitude of the signal component of the 
reconstruction error map is larger for the ML map-making than 
for the destriping. The signal component of the reconstruction 
error depends on A,, in the ML map-making and on Sp in the 
destriping (see Eqs. ( IA.17t and ( IA.19» . The source of both 
quantities is the pixelisation noise p. 

We will first examine the spectrum of the pixelisation noise 
(p = H^p). The PSD of every 60 circle averaged ring of a sim- 
ulated signal-only TOD was calculated and the mean of these 
PSDs was taken over the full one year TOD. The mean PSDs of 
all four simulated signal-only TODs and their associated pix- 
elisation noise streams (p - -Zs, s - TOD) were determined. 
They are shown in Fig. lA.ll 

The level of the pixelisation noise is higher for the TODs 
containing CMB and foreground than for the TODs of CMB 
only. This explains why the signal component of the recon- 
struction error increases when the foreground emissions are in- 
cluded in the simulations. 

The sample pk of the pixelisation noise stream (k indexes 
the sample) is 



Pk 



V * iek / * iek 



(A.21) 



where Sk is the k* sample of the TOD, ; e k refers to those TOD 
samples (including Sk) that hit the same pixel as Sk and A^^ is the 
number of hits in that pixel. We can assume that the magnitudes 
of the pixel-to-pixel correlations are smaller for the tempera- 
ture differences (i,- - Sk) than for the temperatures themselves. 
This means that notable correlation exists between pk and pk- 
only if they are samples from the same pixel. Only the frac- 
tion Npix/Nt I of the TOD samples are observations from 
the same pixel leading to, in average, a weak correlation be- 
tween the samples of the pixelisation noise. This explains the 
"white noise" type nearly flat spectra of the pixelisation noises 
(see Fig. lA.U . The fact that in Fig. lA.lh the pixelisation noise 
is larger for the symmetric beam than for the elliptic beam and 
opposite in Fig. lA.lb reflects the different methods that were 
used in generating the symmetric beam TODs (total convolu- 
tion in case 1 vs. scanning a high resolution map in case 3, see 
Sect.H. 
Vector A„ 



xp is a solution to Eq. (IA.12> . The matrix N is diag- 
onal with samples (bins) of 1 -t- /k// (cf. Eq. (I10» in its diago- 
nals. The diagonal elements of (N - I,) ' are » 1 for frequen- 
cies higher than the knee frequency. Fig. lA.ll suggests that the 
power of H' ZHA^ is considerably smaller than the power of 
A.p leading to an approximation where we can ignore H ' ZH in 
Eq. 



( IA.12> . This indicates that in the ML map-making the sig- 
nal component of the reconstruction error is mainly determined 
by that part of the pixeUsation noise spectrum that falls below 
the knee frequency (Hivon et al.7005'l. 
By looking at Eqs. (0} and iA.l9\ 



we can expect that in 
destriping the uniform baselines of the pixelisation noise con- 
tribute to this error. Because we assumed an exact repetition of 
the pointings of the 60 circles of a ring, the pixelisation noise 
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Frequency [Hz] 



Fig. A.l. The PSD of every 60 circle averaged ring of a simu- 
lated signal-only TOD was calculated (8766 PSDs for a TOD) 
and the mean of the PSDs was taken. The mean PSDs of all 
signal-only TODs generated for this study are shown (curves 
labeled with "TOD"). The mean PSDs of the associated pixeli- 
sation noise streams were calculated as well and are shown too. 
a) TOD contains CMB and foreground, b) TOD contains only 
CMB. The temperature units are antenna yuK at 100 GHz. 



Appendix B: Pixel window and pointing 
distribution effects 

There are a number of effects contributing to the difference be- 
tween the input spectrum C™ and that of the binned noiseless 
map C^. For cases 1, 2 and 4 these include the beam shape, 
spectrum smoothing due to sample interpolation from the total- 
convolver temperature grid (see Sect.|3} and the effective elon- 
gation of the beam due to sample integration. The sample inte- 
gration was simulated in Level S by generating multiple (fast) 
samples at a higher sampling rate and the final output sample 
of the detector (at sampling rate /s) was an average (with equal 
weights) of the fast samples. 

In case 3 the interpolation and sample integration errors do 
not occur, because the signal part of the TOD was made by 
picking the temperatures (at sampling rate /s) from a high- 
resolution input map, which was smoothed by a spherical 
beam. Thus the only beam effect was that of the spherical beam, 
which we already corrected for in Fig. |9l The input spectrum 
displayed in Fig.|9lis actually that of this input map, corrected 
for the beam. 





J 






- C^/(B^cf) 









200 400 600 800 1000 1200 1400 



/ 

Fig. B.l. The ratio of the spectra shown in Fig.|9l(black curve). 
The grey curve shows the same ratio after the spectrum of the 
binned map (C^) has been deconvolved with an approximate 
pixel window of the binning: = D^(512)/D^(1024). 



of those circles is periodic (with period T = 60 s) and it has a 
Fourier series representation whose Fourier mode frequencies 
are multiples of 1 /T. The lowest (zero) mode contributes to the 
uniform baselines, whereas the modes up to the knee frequency 
(0. 1 Hz) contribute to the reconstruction error of the ML map- 
making. This explains why the signal component of the recon- 
struction error is smaller in destriping than in ML map-making. 
The evaluation of the exact effect in the maps is complicated by 
the scanning strategy. 



Another effect comes from how the detector pointings sam- 
ple the sky, or, in this case (case 3), the small pixels of the in- 
put map, to produce the binned map with its larger pixel size 
(the same as the output map). In case 3 the input map had 
A^side = 1024, whereas the binned map had A^side - 512, so 
that each pixel of the binned map can be divided into four sub- 
pixels corresponding to the pixels of the input map. (Note that 
the discussion in Sect.|2on map-making methods assumed the 
same pixel size for input and output maps, and therefore did 
not recognize the effects discussed here.) If each of the 4 sub- 
pixels had been hit by the same number of times, the resulting 
binned noiseless map would be just the input map downgraded 
to A^side =512. The effect on the map spectrum should then be 
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given by the ratio of the two pixel windows D^(512)/D^(1024). 
Here D((512) and Df(1024) are the HEALPix pixel window 
functions fovN^ide = 512 and A^sUe = 1024 (Gorski et al. 2005b i. 
(In the real situation the input map is replaced by the sky with 
"^side = so that the corresponding factor is just the pixel 
window of the binned map.) We show in Fig. lB.ll how this rep- 
resents the effect well up to f ~ 800. 

The remaining effect, which blows up at high i, is due to 
two things: 1) the nonuniform sampling of the four subpix- 
els (or, in the real world, that of the output map pixel area on 
the sky), 2) that the HEALPix pixel window functions them- 
selves represent an approximation, as discussed below. This re- 
maining effect represents coupling between the i modes of the 
spectra, which couples power from the low-^ to the high-/" that 
shows up as an high-^" excess power This effect was discussed 
in Poutanen et al. .2004 where it was modelled as a signal bias. 

Let us examine the distribution of the detector pointings in 
the sky and its impact on the spectrum of the binned noiseless 
map in more detail. The following discussion can be applied 
both to the real case of observing the sky and the case (our case 
3) where the TOD is picked from a high-resolution pixelized 
input map. We consider the samples i, (/ indexes the sample) 
of the CMB-only TOD that fall in a pixel k of the binned (or 
output) map (see Eq. 



iA.lH ). The number of hits in that pixel 
is Nk- We assume that every pixel is hit (100% sky coverage at 
A^side = 512 resolution), so that > 1. The temperature of s, 
can be given as 



Si = ^ atmBcYcniVLi). 



(B.l) 



tm 



Here aim represent the CMB sky (see Sect. 12.31 1. Bt is the re- 
sponse of the symmetric beam and n, is a unit vector pointing in 
the direction of the beam centre (or, in the case where the TOD 
is just picked from an input map, the direction to the centre of 
the input map pixel the detector is pointing at). The temperature 
of the pixel k of the binned map is 



* iek tm * iek 



(B.2) 



where / e k refers to those TOD samples that hit the pixel k. 

The expansion coefficients of the binned map are obtained 
by an inverse spherical harmonic transformation 



(B.3) 



k=0 



We assume a HEALPix pixelisation, where the pixels have the 
same area Qp = An/Npix- The unit vector pointing to the centre 
of the pixel ^ is q^:. 

After inserting from Eq. iB.2\ to Eq. iB.3i we obtain 
for the expansion coefficients of the binned map 

«L = 2 ^e'm'B'fnp 2 F 2 Y{',n'{ndY*Jqk). (B.4) 

t'm' k=0 * i£k 

This equation defines a coupUng matrix 

A'oix-l 



k=0 i.k 



(B.5) 



between the ac,„ of the binned map and the CMB sky. 

Using the statistical isotropy of the CMB sky ((a^mflf/,,,/) = 
6ce'Smm'{Cf)) we obtain for the ensemble mean of the angular 
spectrum of the binned map 

1 ^ 

<C^) = 27TT Z <l«^'"l'> = Z <'^?'<Cr>' (B.6) 

m=-t {' 

where M^^, is the mode coupUng matrix (kernel matrix) of the 
binned map 



Ml = 



1 



2e + \ 



Cmt' 



(B.7) 



In spite of the fact that the binned noiseless map has a 
full sky coverage, its mode coupling matrix is not diag- 
onal but it is only diagonally dominant with small non-zero 
off-diagonal elements, because the pixel area has been nonuni- 
formly sampled (in case 3, hits are in 4 subpixel centres only 
and unevenly distributed among them). The off-diagonal ele- 
ments are responsible for the coupling of the power from the 
low-^ to high-^ that shows up as a high-^ excess power in 
(see Figs. I^andlgm. 

Finally, let us consider what happens if the number of hits 
in the pixel increases, the hits in the pixel area become evenly 
distributed and a symmetric circular pixel shape is assumed. In 
that case the sum ^ Yii£k ^rm'(ii,) can be approximated by an 
integral whose value can be given in a simple form 



y, YcnAni) De,Y(,,„,iqk), 
Nk ^ 



(B.8) 



where Dr ^ Df'(512). (To be precise, the above limiting value 
will be reached, with a different D(, whenever the distribution 
of the hits is the same in every observed pixel and the distribu- 
tion is fully symmetric around its centre q^). Under these as- 
sumptions we obtain an approximation for the mode coupling 
matrix of the binned map 



DiMw. 



(B.9) 



Here Mw is the mode coupling matrix of the MASTER method 
(Hivon et al. l2002l and Sect. l2.3l of this paper). We can see that 
the MASTER approach for the power spectrum estimation (Eq. 
(Q) corresponds to an approximation that a large number of 
hits is symmetrically distributed in every observed pixel. For 
the full sky map (like our binned noiseless map at A^side =512 
resolution) the MASTER mode coupling matrix Mw is close 
to a unit matrix and it cannot explain the high-^ excess power 
that we see in the spectrum of the binned noiseless map. (The 
full sky Mw does have tiny non-zero off-diagonal elements, 
because the spherical harmonics are not exactly an orthogonal 
set of functions in the pixelised sky). 
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